set more 1
clear
set obs 10000

*Q is additional data (not in Clark West)
local Q=0
local R=240
local P=720
local S=`Q'+`R'+`P'



gen time=_n
tsset time

local tsamplestart=`Q'+1
local tfirstpred=`Q'+`R'+1
local tlastpred=`Q'+`R'+`P'

*Define parameters

local vare=18
local varnu=0.025
local r=0.95
local rho=-0.75
local gamma2=0.35

*Define Cholesky Matrix and draw random variables
local c11=`vare'^0.5
local c21=`rho'*(`vare'*`varnu')^0.5/`c11'
local c22=(`varnu'-(`c21'*`c21'))^0.5
di `c11'
di `c21'
di `c22'

forvalues c=1/2{
	gen oosmodel`c'=.
	gen cwmodel`c'=.
	gen cwchange`c'=.
	gen oosr2model`c'=.
}


gen lnr=.
gen dsecond=time<=`Q'+`R'+`P' & time>`Q'+`R'+`P'/2

forvalues j=1/10000{
di `j'

capture drop u1-y2
capture drop pred1-fullmean2
gen u1=invnorm(uniform())
gen u2=invnorm(uniform())

gen e=`c11'*u1
gen nu=`c21'*u1+`c22'*u2


gen y1=0.5+e in 1/`S'
gen z=0.15+(`varnu'/(1-`r'^2))^0.5*nu in 1
replace z=0.15+`r'*z[_n-1]+nu in 2/`S'
gen y2=0.5+`gamma2'*z[_n-1]+e in 1/`S'

gen pred1=.
gen pred2=.
gen modelmean1=.
gen fullmean1=.
gen modelmean2=.
gen fullmean2=.

forvalues i=`tfirstpred'/`tlastpred'{


		local i1=`i'-1
		qui{
		capture drop x
		reg y1 l.z in `tsamplestart'/`i1'
		predict x
		replace pred1=x in `i'
		
		capture drop x
		reg y2 l.z in `tsamplestart'/`i1'
		predict x
		replace pred2=x in `i'
		
		sum y1 in `tsamplestart'/`i1'
		replace modelmean1=r(mean) in `i'
		
				
		sum y2 in `tsamplestart'/`i1'
		replace modelmean2=r(mean) in `i'
		
		
		
		}
		*di `i'
		
}










*R2 Analysis


*returns used below

forvalues c=1/2{

	qui{
	replace lnr=y`c'

	local z=1

		local means `" "model"  "'

		local z1=`z'+1
		local z2=`z'+2
		local z3=`z'+3
		local z4=`z'+4
		local z5=`z'+5
		
		*R2
		capture drop x*
		gen x=(pred`c'-lnr)^2 
		sum x if time>=`tfirstpred' & time<=`tlastpred'
		local sqpred=r(sum)
		
		
		
		foreach mean of local means{
			capture drop x*
			gen x=(`mean'mean`c'-lnr)^2 
			sum x if time>=`tfirstpred' & time<=`tlastpred'
			local sq`mean'=r(sum)
			replace oosr2`mean'`c'=(1-`sqpred'/`sq`mean'')*100 in `z'			
			
		}
		
		
		
		*Tests	
		capture drop x*
		gen xf=(lnr-modelmean`c')^2-((lnr-pred`c')^2-(modelmean`c'-pred`c')^2)				
		reg xf if time>=`tfirstpred' & time<=`tlastpred'
		replace oosr2model`c'=(_coef[_cons]/_se[_cons]) in `z1'
		
		reg xf dsecond if time>=`tfirstpred' & time<=`tlastpred'
		replace oosr2model`c'=(_coef[dsecond]/_se[dsecond]) in `z2'
		
		
		
	
		}

		

	

replace oosmodel`c'=oosr2model`c'[1] in `j'
replace cwmodel`c'=oosr2model`c'[2] in `j'
replace cwchange`c'=oosr2model`c'[3] in `j'
}
}

save sim, replace




*Analyze results
capture drop what size

count if cwmodel1!=.
local N=r(N)
count if cwmodel1>=1.282 & cwmodel1!=.
local x=r(N)
gen what= "size for standard test simulated as in Clark West" in 1
gen size=`x'/`N' in 1


qui count if cwchange1!=.
local N=r(N)
qui count if cwchange1<=-1.282 & cwchange1!=.
local x=r(N)

replace what= "size for second half having lower error" in 2
replace size=`x'/`N' in 2

qui count if cwchange1>=1.282 & cwchange1!=.
local x=r(N)

replace what= "size for second half having larger error" in 3
replace size=`x'/`N' in 3


br what size in 1/3


